import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.signal import find_peaks
from scipy.signal import savgol_filter


# Datos 1
df = pd.read_csv("vext195vret092.csv")
df1 = pd.read_csv("vext195vret192.csv")
df2 = pd.read_csv("vext195vret291.csv")

df["vac"] = pd.to_numeric(df["vac"])
df["va"] = pd.to_numeric(df["va"])
df1["vac"] = pd.to_numeric(df1["vac"])
df1["va"] = pd.to_numeric(df1["va"])
df2["vac"] = pd.to_numeric(df2["vac"])
df2["va"] = pd.to_numeric(df2["va"])

# Visualización de gráficas ---------------------------------------------------

#grafica 1
plt.figure()
plt.locator_params(axis='x', nbins=10)
#plt.scatter(df["vac"], df["va"], label="vret 0.92V")
#plt.scatter(df1["vac"], df1["va"], label="vret 1.92V")
#plt.scatter(df2["vac"], df2["va"], label="vret 2.91V")
plt.plot(df["vac"], df["va"],".", label="vret 0.92V")
plt.plot(df1["vac"], df1["va"],".", label="vret 1.92V")
plt.plot(df2["vac"], df2["va"],".", label="vret 2.91V")

plt.xlabel("Vac (V)")
plt.ylabel("Va (V)")
plt.title("Vext = 1.95 V")
plt.legend()
plt.grid()

# Datos 2
df3 = pd.read_csv("vret150vext150.csv")
df4 = pd.read_csv("vret150vext250.csv")

df3["vac"] = pd.to_numeric(df3["vac"])
df3["va"] = pd.to_numeric(df3["va"])
df4["vac"] = pd.to_numeric(df4["vac"])
df4["va"] = pd.to_numeric(df4["va"])

#grafica 2
plt.figure()
plt.locator_params(axis='x', nbins=10)
#plt.scatter(df3["vac"], df3["va"], label="vext 1.50V")
#plt.scatter(df1["vac"], df1["va"], label="vext 1.95V")
#plt.scatter(df4["vac"], df4["va"], label="vext 2.50V")
plt.plot(df3["vac"], df3["va"],".", label="vext 1.50V")
plt.plot(df1["vac"], df1["va"],".", label="vext 1.95V")
plt.plot(df4["vac"], df4["va"],".", label="vext 2.50V")
plt.xlabel("Vac (V)")
plt.ylabel("Va (V)")
plt.title("Vret = 1.50 V")
plt.legend()
plt.grid()

plt.show()

# Diferencias de V ------------------------------------------------------------

def encontrar_picos(df, nombre_serie, vext, distancia, prominencia):
    indices_max, _ = find_peaks(df["va"], distance=distancia, prominence=prominencia)
    indices_min, _ = find_peaks(-df["va"], distance=distancia, prominence=prominencia)  # <-- nuevo
    
    v_maximos = df["vac"].iloc[indices_max].values
    v_minimos = df["vac"].iloc[indices_min].values  # <-- nuevo
    
    print(f"--- Análisis para {nombre_serie} y vext {vext}V ---")
    print(f"Máximos encontrados en Vac: {v_maximos}")
    print(f"Mínimos encontrados en Vac: {v_minimos}")  # <-- nuevo
    
    DV = []
    for i in range(len(v_maximos)-1):
        DV.append(v_maximos[i+1]-v_maximos[i])
    delta_e_medio = np.mean(DV)
    print(f"Energía de excitación estimada (máximos): {delta_e_medio:.2f} eV\n")
    
    return v_maximos, v_minimos  # <-- modificado

# Ejemplo de uso con uno de tus DataFrames
distancia, prominencia = 10, 0.5
max_df = encontrar_picos(df, "vret 0.92V", 1.95, distancia, prominencia)
max_df1 = encontrar_picos(df1, "vret 1.92V", 1.95, distancia, prominencia)
max_df2 = encontrar_picos(df2, "vret 2.91V", 1.95, distancia, prominencia)

max_df3 = encontrar_picos(df3, "vret 1.50V", 1.50, distancia, 0.12) # Da: 17.80
print("-> da 17.80V")
max_df4 = encontrar_picos(df4, "vret 1.50V", 2.50, distancia, prominencia)


# --- Diferencias de V ajuste parabola ----------------------------------------


def ajuste_parabola_pro(x, y, idx, ventana):
    # Selección dinámica de la ventana
    i_min = int(max(0, idx - ventana))
    i_max = int(min(len(x), idx + ventana + 1))
    
    x_fit, y_fit = x[i_min:i_max], y[i_min:i_max]
    if len(x_fit) < 3: return x[idx], y[idx], None, x_fit, y_fit

    coef = np.polyfit(x_fit, y_fit, 2)
    # Vértice: -b / 2a
    x_ext = -coef[1] / (2 * coef[0])
    y_ext = np.polyval(coef, x_ext)
    
    return x_ext, y_ext, coef, x_fit, y_fit

def analizar_flexible(df, nombre, dist=30, prom=0.01, vent=5, suavizado=11):
    x = df["vac"].values
    y = df["va"].values
    y_suave = savgol_filter(y, suavizado, 3)
    
    indices_max, _ = find_peaks(y_suave, distance=dist, prominence=prom)
    indices_min, _ = find_peaks(-y_suave, distance=dist, prominence=prom)  # <-- nuevo
    
    if len(indices_max) > 3:
        indices_max = sorted(indices_max, key=lambda i: y_suave[i], reverse=True)[:3]
        indices_max = np.sort(indices_max)
    if len(indices_min) > 3:  # <-- nuevo bloque
        indices_min = sorted(indices_min, key=lambda i: y_suave[i])[:3]
        indices_min = np.sort(indices_min)

    v_max, v_min = [], []
    plt.figure(figsize=(8, 4))
    plt.plot(x, y, '.', color='silver', label='Datos')
    
    for idx in indices_max:
        x_e, y_e, coef, x_f, y_f = ajuste_parabola_pro(x, y, idx, vent)
        v_max.append(x_e)
        if coef is not None and coef[0] < 0:
            xf = np.linspace(x_f[0], x_f[-1], 50)
            plt.plot(xf, np.polyval(coef, xf), 'r-', lw=2)
        plt.plot(x_e, y_e, 'ko', label="Máximo")

    for idx in indices_min:  # <-- nuevo bloque
        x_e, y_e, coef, x_f, y_f = ajuste_parabola_pro(x, y, idx, vent)
        v_min.append(x_e)
        if coef is not None and coef[0] > 0:  # parábola hacia arriba para mínimo
            xf = np.linspace(x_f[0], x_f[-1], 50)
            plt.plot(xf, np.polyval(coef, xf), 'b-', lw=2)
        plt.plot(x_e, y_e, 'b^', label="Mínimo")

    v_max = np.array(v_max)
    v_min = np.array(v_min)
    dv_max = np.diff(v_max)
    dv_min = np.diff(v_min)  # <-- nuevo
    dv_m = np.mean(dv_max) if len(dv_max) > 0 else 0
    dv_m_min = np.mean(dv_min) if len(dv_min) > 0 else 0  # <-- nuevo
    phi = v_max[0] - dv_m if len(v_max) > 0 else 0

   # Al final de analizar_flexible, sustituye el print+show por:
    print(f"\n>> {nombre}")
    print(f"   Máximos: {np.round(v_max, 2)} | ΔV: {dv_m:.2f} | Phi: {phi:.2f}")
    print(f"   Mínimos: {np.round(v_min, 2)} | ΔV: {dv_m_min:.2f}")
    
    plt.title(f"{nombre}\n$\Delta V_{{max}}={dv_m:.2f}V,\ \Delta V_{{min}}={dv_m_min:.2f}V,\ \phi={phi:.2f}V$")
    plt.xlabel("Vac (V)")
    plt.ylabel("Va (V)")
    plt.grid(True, alpha=0.3)
    plt.legend(["Datos", "Ajuste máx.", "Máximo", "Ajuste mín.", "Mínimo"])
    plt.show()
    
    return dv_max, dv_min, phi
"""
analizar_flexible(df, "vret: 0.92V, vext: 1.95V", dist=4, prom=0.01, vent=2)
analizar_flexible(df1, "vret: 1.92V, vext: 1.95V", dist=25, prom=0.05, vent=3)
analizar_flexible(df2, "vret: 2.91V, vext: 1.95V", dist=25, prom=0.05, vent=3)
analizar_flexible(df3, "vret: 1.50V, vext: 1.50V", dist=20, prom=0.005, vent=3) 
analizar_flexible(df4, "vret: 1.50V, vext: 2.50V", dist=10, prom=0.1, vent=4)
"""
resultados = []
resultados.append(analizar_flexible(df,  "vret: 0.92V, vext: 1.95V", dist=4,  prom=0.01,  vent=2))
resultados.append(analizar_flexible(df1, "vret: 1.92V, vext: 1.95V", dist=25, prom=0.05,  vent=3))
resultados.append(analizar_flexible(df2, "vret: 2.91V, vext: 1.95V", dist=25, prom=0.05,  vent=3))
resultados.append(analizar_flexible(df3, "vret: 1.50V, vext: 1.50V", dist=20, prom=0.005, vent=3))
resultados.append(analizar_flexible(df4, "vret: 1.50V, vext: 2.50V", dist=10, prom=0.1,   vent=4))

todos_dv_max = np.concatenate([r[0] for r in resultados if len(r[0]) > 0])
todos_dv_min = np.concatenate([r[1] for r in resultados if len(r[1]) > 0])
todos_phi    = np.array([r[2] for r in resultados if r[2] != 0])

# Media global de todos los ΔV juntos (máximos y mínimos)
todos_dv = np.concatenate([todos_dv_max, todos_dv_min])

print("\n========== RESUMEN GLOBAL ==========")
print(f"ΔV medio (máximos):          {np.mean(todos_dv_max):.3f} ± {np.std(todos_dv_max):.3f} V")
print(f"ΔV medio (mínimos):          {np.mean(todos_dv_min):.3f} ± {np.std(todos_dv_min):.3f} V")
print(f"ΔV medio (max + min):        {np.mean(todos_dv):.3f} ± {np.std(todos_dv):.3f} V")
print(f"Potencial de trabajo medio:  {np.mean(todos_phi):.3f} ± {np.std(todos_phi):.3f} V")